#if !defined(IMGLIB_OPERATIONS_H)
#define IMGLIB_OPERATIONS_H

#include <typelib/types.h>
#include <cassert>

namespace img {

/// the function interpolates raster at a given point;
/// if the point lies outside the body,
/// the nearest pixel on the raster's border is returned (function returns false in that case);
/// the function performs transformation of fetched pixels before interpolating them
template <class I, class R> 
bool bilinear(I interpreter, typename I::result_type & res, const R & raster, float x, float y) 
{
	x -= 0.5f;
	y -= 0.5f;

	size_type ux, uy;
	unsigned int code = 0;

	if (x < 0) {
		code |= 1;
		ux = 0;
	} else if (x >+ raster.width() - 1) {
		code |= 2;
		ux = raster.width() - 1;
	}
	else
		ux = size_type(x);

	if (y < 0) {
		code |= 4;
		uy = 0;
	} else if (y >= raster.height() - 1) {
		code |= 8;
		uy = raster.height() - 1;
	}
	else
		uy = size_type(y);

	float ax, ay;

	switch (code) {
	case 0:
		ax = x - ux;
		ay = y - uy;
		res = ::bilinear(interpreter(raster(ux, uy)), interpreter(raster(ux, uy + 1)), interpreter(raster(ux + 1, uy)), interpreter(raster(ux + 1, uy + 1)), ax, ay);
		return true;
	case 4:
	case 8:
		ax = x - ux;
		res = linear(interpreter(raster(ux, uy)), interpreter(raster(ux + 1, uy)), ax);
		return false;
	case 1:
	case 2:
		ay = y - uy;
		res = linear(interpreter(raster(ux, uy)), interpreter(raster(ux, uy + 1)), ay);
		return false;
	case 5:
	case 6:
	case 9:
	case 10:
		res = interpreter(raster(ux, uy));
		return false;
	default:
		res = typename I::result_type();
		return false;
	}		
}

template <class I, class R> 
inline bool bilinear(I interpreter, typename I::result_type & res, const R & raster, float_point p)
	{ return bilinear(interpreter, res, raster, p.x, p.y); }

/// faster version, that uses precalculated int values of coordinates, and does not require +0.5f fox x & y
template <class I, class R> 
bool bilinear_preint(I interpreter, typename I::result_type & res, const R & raster, float x, float y, int ix, int iy) 
{
	//x -= 0.5f;
	//y -= 0.5f;
	assert (x < INT_MIN || x > INT_MAX || floor(x) == ix);
	assert (y < INT_MIN || y > INT_MAX || floor(y) == iy);

	size_type ux, uy;
	unsigned int code = 0;

	if (ix < 0) {
		code |= 1;
		ux = 0;
	} else if (ix + 1 >= raster.width()) {
		code |= 2;
		ux = raster.width() - 1;
	}
	else
		ux = size_type(ix);

	if (iy < 0) {
		code |= 4;
		uy = 0;
	} else if (iy + 1 >= raster.height()) {
		code |= 8;
		uy = raster.height() - 1;
	}
	else
		uy = size_type(iy);

	float ax, ay;

	switch (code) {
	case 0:
		ax = x - ux;
		ay = y - uy;
		res = ::bilinear(interpreter(raster(ux, uy)), interpreter(raster(ux, uy + 1)), interpreter(raster(ux + 1, uy)), interpreter(raster(ux + 1, uy + 1)), ax, ay);
		return true;
	case 4:
	case 8:
		ax = x - ux;
		res = linear(interpreter(raster(ux, uy)), interpreter(raster(ux + 1, uy)), ax);
		return false;
	case 1:
	case 2:
		ay = y - uy;
		res = linear(interpreter(raster(ux, uy)), interpreter(raster(ux, uy + 1)), ay);
		return false;
	case 5:
	case 6:
	case 9:
	case 10:
		res = interpreter(raster(ux, uy));
		return false;
	default:
		res = typename I::result_type();
		return false;
	}		
}

template <class I, class R> 
inline bool bilinear_preint(I interpreter, typename I::result_type & res, const R & raster, float_point p, int_point i)
	{ return bilinear_preint(interpreter, res, raster, p.x, p.y, i.x, i.y); }

template <typename T>
struct cast_interpreter {
	typedef T result_type;
	template <typename U>
	T operator ()(const U & u) const
		{ return T(u); }
};

/// converts pixels to given type before interpolating them
template <class T, class R> 
inline bool bilinear(T & res, const R & raster, float x, float y) 
	{ return bilinear(cast_interpreter<T>(), res, raster, x, y); }

template <class T, class R> 
inline bool bilinear(T & res, const R & raster, float_point p) 
	{ return bilinear(res, raster, p.x, p.y); }

template <class R> 
inline typename R::elem_type bilinear(const R & raster, float x, float y) {
	R::elem_type res;
	bilinear(res, raster, x, y);
	return res;
}

template <class R> 
inline typename R::elem_type bilinear(const R & raster, float_point p) 
	{ return bilinear(raster, p.x, p.y); }

template <class R> 
inline typename R::real_elem_type interpolate(const R & raster, float x, float y) {
	R::real_elem_type res;
	bilinear(res, raster, x, y);
	return res;
}

template <class R> 
inline typename R::real_elem_type interpolate(const R & raster, float_point p)
	{ return interpolate(raster, p.x, p.y); }

} // img

#endif // IMGLIB_OPERATIONS_H
